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Abstract. The two-site Holstein model represents a first non-trivial paradigm for 
the interaction between an itinerant charge with a quantum oscillator, a very common 
topic in different ambits. Exact results can be achieved both analytically and 
numerically, nevertheless it can be useful to compare them with approximate, semi- 
classical techniques in order to highlight the role of quantum effects. In this paper we 
consider the adiabatic limit in which the oscillator is very much slow than the electron. 
A density matrix approach is introduced for studying the charge dynamics and the 
exact results are compared with two different approximations: a Born-Oppenheimer- 
based Static Approximation for the oscillator (SA) and a Quantum-classical (QC) 
dynamics. 
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1. Introduction 

A common problem in chemical reaction and in solid state systems consists in tunneling 
charges between localized sites. The surrounding crystal or molecular system interacts 
with the moving charge and can hugely affect its dynamics. The competition between 
kinetic energy of the charge and localization effects, due to local coupling with lattice 
vibration (phonons), produces the small polaron, i.e. a charge dressed by a cloud of 
multiphonon processes, when the latter prevails [U [2] . 

The two-site cluster is the minimal system in which competition of hopping between 
two sites and phonon localization effects takes place and it has been extensively studied 
since the pioneering work of Holstein |2j. Ground state properties [3j HJ El E] along 
with spectral properties [7J [3j [HJ [H [10] and correlation functions dynamics [HJ [11] have 
been studied using both numerical and analytical methods. Also for its simplicity, 
the two-site cluster has been studied in more involved problems such as polaron 
formation in the presence of double exchange [12J or in the presence of on-site electronic 
repulsion [13]. Recently, this system has also been adopted to explain [TJ] conductance 
experiments in nanotubes [15] The two-site system, introduced in the aforementioned 
mentioned cases of solid state physics, is relevant also in the physics of organic materials 

[ISl EH EH EH [201 [21]. 

In this paper, we study the electron and polaron reduced density matrix in order to 
characterize the transfer dynamics, its degree of coherence and how the temperature 
affects these dynamical properties. In a previous work [9] we used an analytical 
technique, consisting in a mapping into an effective anharmonic oscillator, for studying 
the spectral functions. The same technique can be used also for getting exact results 
for the dynamics. Here, we focus on the adiabatic regime, in which the lattice or 
nuclei dynamics is much slower than that of the electron. In this limit, several 
approximate techniques are usually adopted. Here we compare exact results with a 
Born-Oppenheimer static approach (SA) and a Quantum-Classical (QC) dynamics. The 
comparison between different techniques can be done easily in the two-site system but it 
can be useful also for testing the validity of these methods to a more extended system. 

The model is described in Sect. [2j In Sect. [3] we introduce the reduced density 
matrix for the polaron and the electron. In Sect. H] we describe the exact mapping, 
by means of Fulton-Gouterman transformations, from the original electron-phonon 
problem into two single anharmonic oscillators. Then the QC and SA approximation 
are described. In [5] some result are presented and the comparison between these three 
different techniques is discussed. 

2. The model 

The model describes an electron, in the tight binding approximation, moving in a 
two-site lattice and interacting with it by the local distortion of the lattice site. The 
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Hamiltonian is [5] 

H = —J(c\c 2 +4 c i) +^o{a\ai + a} 2 a 2 ) - g[c\c 1 (a\ + a x ) +4 c 2( a 2 + °2)](l) 

cj and aj are, respectively, the electron and phonon creation operators. The strength of 
the electron-phonon interaction is given by the constant g and J is the electron hopping 
integral. Here we consider one local oscillator per site with frequency ujq. 

We can reduce the the degrees of freedom introducing the coordinates corresponding 
to the center of mass and the relative displacement. The center of mass Hamiltonian 
consists in a displaced oscillator and can be separated from the part depending on 
the relative coordinate (a = (aj — a 2 )/\/2). In the following discussion we shall limit 
ourselves to study only the latter 

H = LUoa^a — Ja x — ga z {a + a). (2) 

In PJ) g = g/ v2 and a pseudo-spin notation has been used by introducing the Pauli 
matrices a z = c\c\ — c 2 c 2 and o x = c\c 2 + c\c\. 

The Hamiltonian (j2J) has a very general form and, even if is has been derived for a 
two-site cluster, it is suitable to describe a very wide class of problems (for example a 
two level system interacting with a single optical mode [7j). 

Beside the temperature, we can choose two parameters that characterize the model 
i) the bare e-ph coupling constant A = g 2 / (uqJ) given by the ratio of the polaron energy 
(E p = —g 2 /uo) to the hopping J and ii) the adiabatic ratio 7 = Uq/J. In terms of these 
parameters we can define a weak-coupling A < 1 and strong coupling A > 1 regimes as 
well as an adiabatic 7 < 1 or anti-adiabatic 7 > 1 regimes. Notice that, in the so called 
atomic (J = 0) limit, the coupling's strength is better described by another constant 
i.e. a = ^/A/(2 7 ) 

In the atomic limit the Hamiltonian is diagonalized by the so-called Lang-Firsov 
(LF) transformation [22] 

D = e aaz{a ^ a) . (3) 

This transformation shifts the phonon operators by a quantity a, while the electron 
operator is transformed into a new fermionic one associated to a quasi-particle, called 
polaron, with energy E p . It can be shown that a 2 is the mean number of phonons in 
the polaron cloud. 

By applying the LF transformation Hq = D^HqD, the atomic Hamiltonian Ho = 

u>oa a — go~ z (a^ + a) becomes 

H = uooa)a + Ep/2, (4) 

the eigenvalues E n = u Q n + E p /2 correspond to the two-fold degenerate eigenvectors 
= D\n,j) = cAri) were the index n = 0, . . . , 00 refers to the photon number, 
j = 1,2 to the electron site number and c! is the polaron creation operator cj = DcjD^ = 
cijexp{(— l) J a(at — a)}. 

In the case of finite J, the hopping term is not diagonalized by ([3]) and the new 
Hamiltonian H = D^HD becomes 

H = uJoa/a — J(o~ x cosh(2a(a^ — a)) + ia y sinh(2a(a^ — a))) + E p /2. (5) 
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Depending on the choice of the parameters, the problem could be better described 
by a electron or polaron excitation picture. In particular, in the weak coupling limit, 
both the small polaron and the electron are good quasiparticles while, in the intermediate 
and strong coupling regimes, the polaron behaviour prevails [23l 19]. 



3. Reduced Density matrix 

Hereafter, we shall assume that charge and oscillator are initially separate, being 
the former localized on the first site and the latter in a mixed thermal state. The 
corresponding density matrix is 

-/3ui n 

p(°) = E— 5H&»W»I®|1><1|, (6) 

n L 

where we used the notation |1) = c\\0). The state \4> n ) depends on the choice of the 
initial preparation [24], in this paper we study two different situations obtained from 
two different limiting regimes: 

(i) electronic preparation (el): The first case corresponds to an initial free Hamiltonian 
(g = 0) where the oscillator is at its thermal equilibrium. In this case \4> n ) = |n). 

(ii) polaronic preparation (pol): In the second case we start from the atomic limit 
(J = 0). Here the oscillator is displaced and the basis in which it thermalizes is 
\4>n) = \ipn) i- e - ^ ne displaced oscillator eigenbasis. 

The dynamics is obtained by switching on g, in the first case, and J, in the second 
one, and letting evolve the density matrix with the Hamiltonian (T2]) p(t) = e~ lHt p(0)e lHt . 

Tracing over the oscillator degree of freedom, we obtain the electron reduced density 
matrix 

p«(t) = Tr ph {p(t)}, (7) 
which, in terms of the oscillator's number states, is 

-/3iv n 

P {d) (t) = J2^^(m\e- im \n,l)(n,l\e iHt \m), (8) 

n,m 

To characterize the motion of the polaron we cannot reduce the density matrix by 
tracing out the phonon degrees of freedom, this is because the polaron itself contains 
phonons. In order to understand better the polaron dynamics, let us first apply a 
Lang-Firsov transformation, the new fermionic particle corresponds to a polaron, so the 
density matrix with the initial localized polaron can be written as 

p<*° l \t) = Tr ph {Dip{t)D R }, (9) 

and reads in terms of the oscillator's number state as 

-f3u) n 

pipol)^ = £ ^_( m | e - im |n, l)(n, l\e lHt \m). (10) 

n,m 

The diagonal elements of the reduced density matrix, in the site basis, represent 
the population of each site while the off-diagonal elements are called coherences 
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and represent the quantum interference between localized amplitudes. An important 
quantity, which is also base-independent, is the so-called purity defined as 

P = Tr{p 2 }. (11) 

It is easy to see that 1/2 < P < 1 and it measures how much the state is pure (P = 1 
if and only if the state is pure and P = 1/2 when it is maximally mixed). 



4. Analytical approaches 

4-1. The adiabatic, static, approximation 

The case, in which a light quantum particle interacts with much more massive particles, 
is very common in solid state and molecular physics. We discuss the adiabatic regime, 
meaning that in a characteristic time for the the light particle dynamics the heavy 
degrees of freedom can be considered approximately quiet. Here we describe the SA 
approach in its basic formulation for the dynamics. 

The Hamiltonian (J2]) can be written in the coordinate-momentum representation 

p 2 mu 2 2 g 

H — x 7=xa z — Ja x — w , (12) 

2m 2 ^ v ; 



with g = g^/lmujQ. In the adiabatic limit (7 1), the phonon is much slower than the 
electron (heavy phonon and large electron tunnelling amplitude) and one can neglect the 
phonon kinetic term in ffT2l . This is the well known Born-Oppenheimer approximation. 
In practice, it consists in studying the electronic problem with x as a classical parameter. 



Within this approximation, we put ujq = yk/m — > (m is the ion mass) and the 
Hamiltonian becomes 

Had = t^x 2 ~ -^ X(Tz ~ Ja *- ( 13 ) 

The eigenvalues can be expressed trough the classical displacement x 

V±(x) = ^x 2 ±Q(x), (14) 



with Q(x) = \j\x 2 + J 2 . The lowest branch (— ) of flHJ) defines an adiabatic potential 
which has a minimum at x = as far as A < 1 while for A > 1, it becomes double 
well potential with minima at x = ±x m = ±\J — 2 4^- 1 in this case the electron is 
mostly localized on a given site. The quantum fluctuations are able to restore the 
symmetry in analogy to what happens for an infinite lattice [25]. It is worth noticing 
that, in this limit, Hamiltonian (J2J) is equivalent to the adiabatic version of the spin- 
boson Hamiltonian [26, 27J. 

The temporal evolution is given by 

e -iH ad t = e -i h £t[ C0S fi( x )t + i( J* a z + -l—cr x ) smQ(x)t], (15) 

\/2U[x) \l[x) 

so the density matrix dynamics can be explicitly calculated 
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The electronic initial preparation, corresponds to the density matrix 

p(0) = |l)(l|^|cia;e-^ 2 |x)^l, (16) 
tracing out the phonon we obtain the electron reduced density matrix with elements 



(el) 
P2,2 



2tt J """" 1 + A 2 M 2 



= -n/^F / d«e-^^ ^_ ^ (17) 



2tt J 2VX 



2 U 2 



where the scaled lenght u = xk^/2/g was introduced. 

In the same way we can introduce the polaronic preparation 

p(0) = |l)(l| e -^/M / dxe-^ x2 -^ x) \x)(x\, (18) 
V 2n J 

It is worth stressing that, in the adiabatic limit, we cannot define the polaronic dynamics, 
as introduced in ([9]), because the operator D is not defined for u = 0, so we better 
have to talk about an electronic dynamics with an initial polaronic preparation. The 
corresponding reduced density matrix is 



Apol) 

P%2 



Jpot) 
Pl,2 



l pj\ 

J- 



due /3;2 *• 



I/3JX 



du e 



u _ 1)2 s m 2 (Jt^(u\) 2 + l) 
(uX) 2 + 1 
uX sin 2 (JtJ{Xu) 2 



((Xu) 2 + 1) 



— I- 



sm(2Jt^J(Xu) 2 + 1) 
2 v /(Am) 2 + 1 



(19) 



It is possible to demonstrate that p'f ^ is actually the adiabatic limit of the diagonal 
element of the reduced polaronic density matrix, while this is not true for the off-diagonal 
elements. 

We want stress that, in the SA approach, the phonon is completely static because 
its momentum p has been neglected. Here, only the initial phonon distribution plays a 
role, but during electron hopping, oscillator is taken to be fixed. 



4-2. A quantum- classical dynamics approximation 

In the adiabatic limit, the slow variables can be considered as classical and a mixed 
quantum-classical dynamics can be introduced. In the past, a lot of schemes for 
quantum-classical dynamics has been proposed, for example starting from the Born- 
Oppenheimer (SA) adiabatic approximation for the ground state at each step and using 
a density functional Hamiltonian [281 122] • Another approach, good for a short time 
dynamics, consists in a mapping from the Heisenberg equations to a classical evolution 
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by an average over the initial condition [30, 3 lj . Some schemes are based on the evolution 
of the density matrix coupled to a classical bath [32], [33]. A systematic expansion over the 
mass ratio has also been done, starting from partial Wigner transform of the Liouville 
operator, in [3U [351 EH]- The QC approximation we use is essentially that of refs. 
[321 [33]. 

Let us now consider Hamiltonian ffl2|) , where x and p are now two classical variables 
which can be represented as the components of a vector u. Then a QC vector can be 
introduced 

/ x \ 

P 



U (X) 0" 



0~n 



(20) 



V <?z J 

The classical variables evolve with the Ehrenfest equations, defined for the averages 
quantities 



in 



x 

while the quantum variables evolves in the Heisenberg picture 
& x = -\f2gxOy 



(21) 



& y = \[2gxo x — 2 Ja z . 



(22) 



Ua,„ 



To give a unified description of the overall evolution, we define a Liouvillian operator 

/ -V2gx \ 



£ = L x + C p + C-a with 



V2gx -2J 
\ 2J J 



(23) 



and C x = x-jr- L v = p4-. So, the time evolution is given by 



' dx P " dp ' 
iCt, 



v(t) = e l "v(0). (24) 

The numerical integration can be implemented using the symetrized Trotter 
breakup formula [3TI [38] 

v(t) ~ (e i£ ^e i£ ^e i£ ^e i£ ^e i£ ^) Af v(0) (25) 

with e = t/N. All the density matrix elements, can be expressed in terms of elements 
of v(t). 

This approach, in contrast with SA, takes into account the dynamics of the phonon, 
even if it remains classical. 
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4-3. Exact diagonalization 

As shown by Fulton and Gouterman |39j, a two-level system coupled to an oscillatory 
system in such a manner that the total Hamiltonian displays a reflection symmetry, may 
be subjected to a unitary transformation which diagonalizes the system with respect to 
the two-level subsystem [lOl HT1 139] . This method can be generalized to the N-site 
situation, if the symmetry of the system is governed by an Abelian group [¥T] . 

In particular, an analytic method for calculating the Green functions of the two-site 
Holstein model is given in [T2J [9]. Here the Hamiltonian is diagonalized in the fermion 
subspace by applying a Fulton Gouterman (FG) transformation. So the initial problem 
is mapped into an effective anharmonic oscillator model. It is possible to introduce 
different FG transformations for the electron and the polaron. The resulting problem 
results to be very simplified and very suitable to be numerically implemented. Analytical 
continued-fraction results, exist for the electron case [91 [12]. 

In this section, we briefly recall the FG transformations method. The density 
matrix elements are given explicitly in terms of effective Hamiltonians and calculated 
by means of exact diagonalization 

The FG transformation we use for the electronic case is 



the new Hamiltonian H = VHV 1 becomes diagonal in the electron subspace 



the diagonal elements, corresponding to the bonding and antibonding sectors of the 
electron subspace, being two purely phononic Hamiltonians 



The operator (—1)" a is the reflection operator in the vibrational subspace and it satisfies 
the condition (— l) ata a(— l) ata = —a. A wide study of the eigenvalue problem was 
carried out in [12] both numerically and analytically by a variational method, extending 
the former results given in |43j. In [12] H± is approximately diagonalized by applying a 
displacement, the dynamics is reconstructed by the calculated eigenvectors and energies. 

The evaluation of the polaron Green function can be done on the same footings 
but the expression involves also the non diagonal elements of the resolvent operators 
causing an exponential increasing of the numerical calculations. 

To avoid this problem, we first perform the FG transformation and then apply, on 
the resulting Hamiltonian ([5]), a different FG transformation 




(26) 




(27) 



H± = u a *a =)= «/(— 1) 



g(a j + a). 



(28) 




(29) 



The new Hamiltonian H 



LF 



V,HV,-' is 




(30) 
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is real and symmetric but not tridiagonal in the basis of the harmonic oscillator, the 
matrix elements of H± are given in [9]. 

In order to write down the density matrix elements, let us introduce the following 
notation: 



R£l(t) = <m|e-^|n>, 

JV# n (t) = (m,l\e- iHt \n,l) = 
= (-l) n + m N™f(t) 

N™f(t) = (m,2\e- iHt \n,l) = 
= (-l) n+m N^(t), 



2 

(-1) 



(32) 
(33) 



(34) 



(35) 



R^nit) + (-l) m+n ^i(i)l (36) 



-iHt 



IVi 1) = j [(-l) n ^L~iW - (-l)^SW] • (37) 



The reduced electron density matrix elements are 



\N^ n (t)\' 



(38) 



the calculation for the polaron case gives 



= E 



in *-<m. n / . \ i 2 



|MT\' ft (*)|' 



fti?(t) = E ^M 2 y(t)M*r(t). 



-/3aJon 



(39) 



5. Comparison between different techniques 

In this section, we show a comparison between the results obtained in the three different 
ways described before: exact diagonalization (ED) by means of mapping introduced in 
Sect. 14.31 the quantum-classical (QC) dynamics approach described in Sect. 14.21 and 
the static Born-Oppenheimer (SA) approximation (Sect. 14.11) . We shall limit ourselves 
to an adiabatic case (7 = 0.1) with electron-phonon interaction strong enough to allow 
the polaron formation (A = 2). 

Let us start from the electron case (left panels in Fig. (pQ)). At short times, 
we see that the state, not being an an eigenstate of H, starts to oscillate coherently 
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Figure 1. Electron (left panels) and Polaron (right panels) dynamics at different 
temperatures T/ujq = 0.1; f; 10. Adiabatic and strong coupling regime^ = 0.1 and 
A = 2. For each temperature are shown the transition probability and the purity. 



with frequency J, as expected for a not interacting electron. Because of interaction, 
oscillations are damped and the electron starts to be entangled with phonon. The time- 
range considered here is of the order of Uq 1 . On this timescale we can see that the higher 
the temperature, the higher the damping of the transition probability and, on the same 
time, the coherence goes down. This means that the electronic state tends to localize 
on one of the two sites with an equal classical probability, in other words it becomes 
a classical mixture of localized states. This is an interesting example of system going 
to a sort of equilibrium despite its finite size. This is because, in the adiabatic limit, 
the oscilltor spectrum tends to a continuum. Actually, since the system is not exactly 
adiabatic, a partial recoherence is gained over much larger timescales of the order of the 
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Figure 2. Time average transition probability as function of time in the SA 
approximation. 

level spacing of the interactting system. Analysis of this regime is beyond the scope of 
the present paper and will appear in a forthcoming publication. Let us now compare the 
three different techniques. At low temperature the classical phonon is almost freezed, 
and so both SA and QC approximation are equivalent. Nevertheless, the ED behaviour 
is quite different because in the exact dynamics also quantum fluctuations play a role. 
A transition to a semiclassical behaviour can be seen for T is comparable with u when 
SA and QC reproduce the ED dynamics. Notice that the same occurs in the transport 
of extended system where classical incoherent transport is achieved when T is greater 
than 0.2lj [H] so that when T ~ wo we have a classical incoherent transport. 

At high temperature the oscillator dynamics plays a relevant role, QC is a much 
better approximation of ED than SA. This fact can be understood by noting that the 
main temperature effect is the damping of the coherent tunnelling oscillations. Once 
these oscillation are sufficiently suppressed, the phonon driven dynamics prevails. In 
the SA framework, the initial thermal distribution of the phonon coordinate makes 
the electron thermalizes irreversibly, in a time that is the shorter the greater the 
temperature. Before this adiabatic thermalization, i.e, in a tunnelling period, the SA 
is still a good approximation. Afterwards, the tunnelling dynamics saturates and the 
phonon dynamics has to be taken into account. 

As far as the polaron is concerned, we can see (right panels in Fig. ([!])), that, 
for low temperature, the polaron is still a good picture, keeping the state quite pure 
even at long times. Nevertheless, the polaron is very localized and its transition 
probability is extremely low. The polaron is trapped inside the initial site. As the 
temperature increases, state becomes mixed but, in contrast with the electron, the 
transition probability increases. So, in the adiabatic limit, the polaron mobility is 
increased by a thermal activation. That can be also seen in Fig. ([2]), where we reported 
the SA time averages of the probability for the charge to be on site 2, in both electron 
and polaron cases. As one can see, the curve is monotonic decreasing for the electron, 
while the polaron mobility starts to be enhanced by the temperature before going down 
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as for the electron. 

As it was stressed before, SA ad QC do not provide a good limit for the purity, as 
we can verify by looking the right panels in Fig. (CQ)). Nevertheless, for the population 
dynamics, are still valid the same considerations done for the electron. 

6. conclusion 

We introduced a density matrix approach for characterizing the reduced charge dynamics 

in a two-site Holstein model. Both the single site population and the purity are 

studied in order to connect the charge transfer with its coherence. In the adiabatic 

limit, two common approximation has been compared with exact results, the aim is 

to highlight the limits of validity of these approximations and to provide a simple 

testing tool, the two-site model, for generalizations to other extended models. Electron 

and polaron dynamics has been compared. Despite the two-level description of the 

moving charge, temperature induces decoherence on intermediate timescale, due to the 

strong interaction with the oscillator. In the polaron case, temperature can enhance the 

(incoherent) charge transfer. 
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